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Abstract-  In  this  paper  we  propose  a  method  for  joint  estimation 
of  delay  and  scale  factor  between  two  deterministic  and 
unknown  signals.  The  method  is  based  on  the  separation  of  the 
two  parameters.  The  scale  factor  is  estimated  from  the 
autocorrelation  functions  (ACFs)  of  the  two  signals.  The  ACFs 
are  independent  on  the  delay  while  they  maintain  the  scale 
factor.  The  delay  is  estimated  by  scaling  one  of  the  two  signals 
by  the  estimated  scale  factor,  so  that  the  delay  estimation  is  not 
biased  by  the  scale  factor.  The  method  proposed  is  compared 
with  the  maximum  likelihood  joint  estimation  and  provides 
better  performance  for  a  rather  large  range  of  signal  to  noise 
ratios.  Potential  applications  of  the  technique  in  surface  EMG 
signal  analysis  are  discussed  and  results  related  to  the  M-wave 
processing  are  reported. 
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estimators,  motor  unit  conduction  velocity  distribution,  M-wave, 
scale  factor  estimation 

I.  Introduction 

Many  works  are  reported  in  literature  in  the  area  of  Radar 
and  passive  Sonar  to  jointly  estimate  the  scale  and  time  delay 
between  two  or  more  signals.  The  most  common  approach  is 
the  maximization  of  the  ambiguity  function  of  the  two 
parameters  [6].  The  problem  usually  encountered  is  the  fact 
that  the  ambiguity  function  depends  on  signal  type  and  model 
[6] .  Several  authors  used  the  decorrelation  method  to  estimate 
the  two  parameters  but  this  approach  needs  a  long 
observation  window  [9]. 

Surface  EMG  signals  detected  during  electrically  elicited 
muscle  contractions  represent  the  synchronous  summation  of 
the  action  potentials  generated  by  the  active  muscle  fibers. 
These  signals  are  quasi  periodic  since  at  each  electrical 
stimulus  a  compound  potential  (called  M-wave)  is  generated. 
The  scale  factor  and  shape  of  each  M-wave  depend  on  the 
conduction  velocity  (CV)  distribution  of  the  active  Motor 
Units  (MUs)  and  on  the  properties  of  the  intra-cellular  action 
potentials.  The  source  characteristics  may  change  with 
fatigue  and  this  is  usually  monitored  by  spectral  analysis  of 
the  M-waves  or  estimation  of  CV  by  estimating  the  delay 
between  two  signals  in  two  different  locations  along  the 
muscle.  Spectral  variables  and  CV  estimated  during  time  in 
sustained  electrically  elicited  contractions  provide  some  basic 
indication  about  the  decrease  in  velocity  of  propagation  of  the 
muscle  fibers.  Since  in  this  case  the  central  nervous  system 
does  not  regulate  the  firings  of  the  MUs,  this  acquisition 
modality  is  often  used  to  investigate  peripheral  phenomena  of 
muscle  fatigue. 

In  this  work  we  propose  a  novel  approach  to  estimate 
jointly  the  time  scale  factor  and  time  delay.  The  method  is 
based  on  a  first  step,  which  aims  at  the  estimation  of  the  scale 
factor  independently  from  the  time  delay,  and  a  second  step 


in  which  the  scale  factor  between  the  two  signals  is 
compensated  by  dilatation  or  compression  of  one  of  the  two 
signals  and  the  time  delay  is  estimated  without  the  bias  due  to 
the  scale  factor.  Potential  applications  of  this  method  to 
surface  EMG  signal  analysis  are  proposed  and  some  results 
on  M-wave  analysis  are  reported. 

II.  Methodology 

A.  Parameter  separation 

An  analytical  model  of  two  signals  delayed  and  scaled  is 
the  following: 

Xi(t)  =  s(t)  +  Wi(t) 

^2(1)  =  s(ca  -d)  +  W2(t) 

where  a  and  d  are  the  scale  factor  and  time  delay, 
respectively,  and  Wj(t)  and  W2(t)  are  independent,  white,  zero 
mean,  additive  gaussian  noises  with  equal  variance  o^. 

A  solution  for  the  problem  of  joint  estimation  of  scale  and 
delay  factor  is  to  separate  the  two  parameters  (scale  and 
delay),  to  estimate  one  of  them  first  and  then  to  estimate  the 
second  after  correction  for  the  first  parameter.  To  estimate  the 
parameters  separately,  the  first  step  consists  in  the  estimation 
of  the  scale  factor  using  the  temporal  autocorrelation 
functions  (ACEs)  of  the  signals.  Indeed,  the  ACE  t)  of  a 
signal  preserves  the  scale  factor  and  is  insensitive  to  the  time 
delay.  Eor  the  model  in  (1)  we  get: 

" "  a 

where  5  is  the  Dirac  distribution,  is  the  white  noise  power, 
and  is  the  temporal  autocorrelation  function  of  s(t). 

The  effect  of  white  noise  can  be  reduced  by  neglecting  the 
value  of  ^  interpolating  this  point 

using  the  neighbouring  points  and  the  best  second  order 
polynomial  fit.  To  get  ride  of  the  multiplication  factor  1/a,  an 
amplitude  normalization  is  then  carried  out.  After  the 
interpolation  and  normalization  we  get  (in  discrete  time 
domain): 

A  1  ^ 

Rx2X2  (/)  =  Rss  (jyV) 

a 

where  ^  is  for  estimated  values.  In  the  same  way  we  obtain: 

A  A 

Rx^x^  (1)  =  Rss(l) 

To  estimate  the  scale  factor  between  R^  ^  (/)  and  R^  „  (/)  , 

Y  ^  X2X2 

the  value  of  a  leading  to  the  minimum  mean  square  error 
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between  the  two  scaled  ACFs  can  be  used.  A  problem  of  this 
method  is  the  limited  estimation  resolution.  To  obtain  high 
resolution,  small  steps  for  a  must  be  in  fact  chosen. 
Resolution  can  be  enhanced  by  interpolating  the  mean  square 
error  function  around  the  minimum  by  a  second  order 
polynomial  approximation  but  this  interpolation  requires  the 
symmetry  of  the  function  around  its  minimum  and  it  can  be 
shown  that  the  mean  square  error  function  is  not  symmetric  in 
this  case  [8].  To  avoid  the  problems  of  calculation  time  and 
loss  of  resolution,  other  methods  using  the  scale  domain, 
discussed  below,  can  be  proposed. 

B.  Scale  transform 

The  scale  transform  of  a  signal,  introduced  firstly  by  Cohen 
[1],  provides  the  representation  of  a  support  limited  signal  s(t) 
in  the  scale  domain.  The  scale  transform  is  defined  as: 


After  a  change  of  variable  t  — >  Ate"  ,  the  transform  in  Eq. 
(2)  can  be  rewritten  as  [8]: 

D,(c)  =  e‘‘'^s(Ate‘‘)e-^^"du  + 

(3) 

where  the  first  integral  corresponds  to  the  period  of  the  signal 
s(t)  between  s(0)  and  the  first  sample  s(l)  of  the  sampled 
version  of  the  signal  s(t).  In  Eq.  (3)  the  terms  between 
parenthesis  show  that  the  scale  transform  is  the  Eourier 
Transform  of  the  exponentially  resampled  and  weighted 
version  of  the  signal  s(t).  Moreover,  if  the  signal  s(t)  is  scaled 
and  energy  normalized  (see  Eq.  4  below),  the  scale  transform 
is  equivalent  to  the  Eourier  transform  of  a  new  function  S  that 
is  translated  by  the  amount  Ina  (Eq.  5).  This  new  translated 
function  is  obtained  by  the  exponential  resampling  and 
weighting: 

Qxp  .resampling 

4as{ca)  ^  +  lna)  (4) 

The  Eourier  transform  of  Eq.  (4)  is: 

5(f  +  lna)Se'^"“D,(c)  (5) 

where  F.T.  stands  for  the  Eourier  transform. 

C.  Estimation  of  scale  factor  from  the  scale  transforms  of  the 
autocorrelation  functions  of  the  two  signals 

Eq.  (5)  shows  the  possibility  of  estimating  the  scale  factor 
a  by  a  spectral  matching  technique  in  the  scale  domain 
between  a  reference  signal  and  its  scaled  version.  The 
advantage  of  this  method  is  the  high  resolution  estimate  of  a 
and  the  rapid  calculation  time  by  using  gradient  methods  for 
the  spectral  matching  [7].  A  drawback  of  this  method  is  the 
presence  of  local  minima.  To  avoid  local  minima  we  can 
choose  a  starting  value  (coarse  value)  of  a  close  to  the  exact 
value.  The  coarse  estimation  can  be  obtained  by  the 
maximum  point  of  the  cross-correlation  function  between  the 
exponentially  resampled  and  weighted  scale  transforms. 


After  scale  estimation,  the  scale  factor  can  be 
compensated  by  scaling  one  of  the  two  signals,  for  example 
by  inverting  the  scale  transform  or  by  using  sine  functions. 
Then  the  delay  is  estimated  as  the  maximum  point  of  the 
cross-correlation  function  (with  a  parabola  fitting  the  apex). 

Eig.  1  shows  the  schematic  representation  of  the  complete 
method  used  for  the  joint  scale-delay  estimation. 


Fig.  I.  Schematic  representation  of  the  method  for  joint  time  delay  and  scale 
estimation.  See  text  for  details. 


III.  VALIDATION  OE  THE  METHOD  ON  SIMULATED 
SIGNALS  AND  POTENTIAL  APPLICATIONS  IN 
SUREACE  EMG  SIGNAL  ANALYSIS 

A.  Validation  of  the  method  on  synthetic  single  MUAP  signals 

Eig.  2  shows  the  results  of  the  application  of  the  method 
on  second  order  Hermite-Rodriguez  functions  corrupted  by 
additive  white  gaussian  noise  (scale  estimation  (a),  and  time 
delay  estimation  (b)).  The  method  proposed  is  compared  with 
the  exhaustive  Maximum  Likelihood  (ML)  and  the  Cramer- 
Rao  lower  bound  (CRLB). 

B.  Potential  applications  for  surface  EMG  signal  analysis 

The  estimation  of  scale  and  delay  may  have  useful 
applications  in  surface  EMG  signal  analysis  both  during 
voluntary  and  electrically  elicited  contractions.  The  detection 
of  single  MU  multi-channel  EMG  signals  during  voluntary 
contractions  allows  estimation  of  CV  of  individual  MUs 
[2]  [3].  The  estimation  of  the  scale  factor  between  MUAPs 
detected  by  different  electrode  systems  along  the  muscle 
fibers  may  give  an  indication  of  fiber  inclination  with  respect 
to  the  detection  array  and  may  be  useful  for  electrode 
placement  or  for  compensating  the  bias  in  CV  estimation  due 
to  the  inclination  angle.  The  same  technique  can  be  used  for 
precise  electrode  placement  when  a  real  time  averaging  of  the 
surface  detected  EMG  signals  is  performed  with  intra¬ 
muscular  detected  potentials  as  triggers  [5]. 


indication  of  mean  CV  if  the  estimation  is  not  biased  by  the 
scale  factor. 


Fig.  2.  Variance  of  estimation  of  the  scale  factor  and  delay  with  the  proposed 
method  and  the  Maximum  Likelihood  joint  estimation  for  different  signal  to 
noise  ratios.  The  Cramer  Rao  lower  bound  is  also  shown.  The  simulated 
signals  are  Hermite-Rodriguez  functions  of  the  second  order  which 
approximate  real  MUAPs  detected  in  double  differential  mode.  The  MUAP 
duration  is  approximately  12  ms  and  the  sampling  rate  is  2048  Hz.  The  scale 
factor  is  a=0.9  and  the  time  delay  is  of  6  samples. 

In  case  of  electrically  elicited  contractions,  the  estimation 
of  scale  factor  between  consecutive  in  time  M-waves  gives 
indications  about  the  change  of  mean  CV,  as  alternative  to 
spectral  analysis  or  CV  estimation.  Moreover,  since  the  scale 
factor  of  an  M-wave  depends  on  the  CV  distribution,  it  can  be 
used  as  indicative  of  fiber  type.  In  particular,  detecting  the 
EMG  signal  along  the  muscle  fibers  during  stimulation,  the 
different  M-waves  will  have  different  scale  factors  due  to  the 
increasing  spread  of  the  compound  potential  as  a 
consequence  of  the  different  time  delays  of  the  MUAPs 
constituting  the  M-wave.  Thus  estimation  of  the  time  delay 
between  compound  potentials  during  time  gives  indication 
about  mean  CV  and  estimation  of  scale  factor  between  M- 
waves  detected  at  the  same  instant  of  time  in  two  different 
locations  in  space  gives  indication  of  CV  distribution  spread 
and  its  changing  during  time.  Moreover,  the  estimation  of 
delay  between  consecutive  in  time  M-waves  gives  also  an 


C.  Simulation  of  electrically  elicited  surface  EMG  signals 

EMG  signals  detected  during  electrically  elicited 
contractions  were  simulated  with  the  model  proposed  by 
Earina  &  Merletti  [4].  MUs  were  placed  randomly  in  the 
detection  volume.  Simulated  signals  were  detected  with  three 
single  differential  recording  systems  placed  along  the  fiber 
direction.  Two  double  differential  signals  were  computed  as 
the  difference  between  consecutive  single  differential 
recordings.  Different  values  of  CV  distribution  mean  and 
standard  deviation  were  simulated.  The  scale  factor  between 
M-waves  detected  in  double  differential  mode  from 
consecutive  recording  systems  was  estimated.  As  expected,  it 
has  been  found  that  the  scale  factor  estimation  is  insensitive 
from  the  thickness  of  subcutaneous  tissue  layers.  On  the  other 
hand,  the  scale  factor  indicates  the  different  CV  distribution 
standard  deviations.  Statistically  different  values  of  scale 
factors  between  M-waves  have  been  found  for  differences  in 
CV  distribution  standard  deviation  lower  than  0.03  m/s.  Eig. 
3  reports  the  results  obtained  in  one  simulation  set. 


CV  distribution  standard  deviation  (m/s) 

Fig.  3.  Scale  factor  (mean  ±  standard  deviation)  between  simulated  M-waves 
detected  in  single  differential  mode  at  a  distance  of  1  cm  along  the  muscle 
fibers.  Fifty  simulations  have  been  performed  for  each  CV  distribution 
standard  deviation  value.  In  each  simulation  the  location  of  the  MUs  in  the 
muscle  was  randomly  changed  in  the  detection  volume  in  order  to  assess  the 
variability  of  the  estimates  due  to  the  volume  conductor  masking  effect. 

V.  Discussion 

The  estimation  of  delay  and  scale  factor  between  two 
deterministic  signals  is  not  trivial  since  the  two  parameters 
are  correlated  and  a  joint  estimation  is  necessary.  We  have 
proposed  a  method  for  the  joint  estimation  which  allows  to 
reach  good  performance  with  limited  computational  time. 
The  performances  of  the  method  are  better  than  those  of  an 
exhaustive  ML  based  estimation  even  if  we  take  a  high 
resolution  search  grid  on  the  ambiguity  function. 

The  technique  can  be  applied  in  many  fields.  Potential 
applications  to  surface  EMG  signal  analysis  have  been 
discussed  and  results  on  simulated  electrically  elicited 
contractions  have  been  reported.  In  particular,  the  estimation 


of  the  scale  factor  between  M-waves  detected  along  the 
muscle  fibers  may  give  indication  about  the  distribution  of 
CVs  with  which  the  action  potentials  propagate.  Since  the 
measure  is  not  absolute  but  relative  (scale  difference  between 
two  M-waves),  it  is  insensitive  to  the  thickness  of 
subcutaneous  tissue  layers,  to  mean  CV  or  to  the  length  of  the 
intra-cellular  action  potentials;  these  factors,  in  fact,  affect  in 
the  same  way  the  width  of  both  M-waves.  On  the  other  hand 
this  estimate  may  be  affected  by  fiber  inclination  with  respect 
to  the  detection  system. 

Although  the  possibility  of  quantitatively  associate  to  a 
scale  factor  estimate  a  corresponding  CV  distribution 
standard  deviation  value  should  still  be  investigated,  the 
presented  results  shows  the  feasibility  of  using  the  scale 
factor  as  a  variable  for  relative  comparisons  between  subjects 
and  muscles.  A  validation  of  the  technique  on  experimental 
signals  will  be  the  subject  of  further  investigation. 

VI.  Conclusions 

A  method  for  joint  scale  and  delay  estimation  and  its 
potential  applications  to  the  analysis  of  electrically  elicited 
EMG  signals  has  been  proposed.  Performances  have 
been  statistically  analysed  (Monte  Carlo  method)  giving 
better  result  in  comparison  to  ML  approach.  In  order  to 
investigate  the  feasibility  of  applying  the  proposed  method  to 
real  cases,  simulated  signals  have  been  generated  leading  to 
promising  results. 
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